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Introduction 

The problem of generating grids for flow calcu- 
lations over airplane configurations is important be- 
cause the grid represents the object on which the 
flow code operates. Several sophisticated procedures 
have been developed for generating grids for fairly 
complex configurations from specifications read in 
from the screen by interactive methods. (See ref. 1.) 
However, such hands-on interactive methods are too 
time-consuming to be feasible when geometries must 
be modified frequently and automatically, as in opti- 
mal design calculations. Many trade-off studies can 
be performed without the complete detailed airplane 
geometry, as for example, wing-fuselage or wing- 
fuselage-fin configurations. For such simplified con- 
figurations, great generality in grid-generation tech- 
niques is not required and simpler procedures that 
can be automated become appropriate. 

The purpose of this paper is to describe one such 
automatic grid-generation procedure. The algebraic 
method used is a variant and generalization of the 
method of transfinite interpolation. Examples of its 
application to generating an Euler-marching grid and 
a grid for subsonic flow calculation are provided. 

Symbols 

A cross-section aspect ratio 

i,j,k indices for £,77, £ grid surfaces, respectively 

p input parameter to control x 

q superellipticity exponent (eqs. (3) and (4)) 

r radius of grid outer boundary 

r two-dimensional position vector with 

components y, z 

s arc length 

x, y, z Cartesian coordinates (longitudinal, lateral, 
and vertical body axis, respectively) 

z z location of symmetry line for superellipse 

a, 7, 6 blending functions 

6 superellipse angle parameter 

A scaling function 

£, 77, £ grid coordinates 

p superellipse scale parameter 

X superellipse eccentricity parameter 

ip variable defined by equation (8) 


Subscripts: 

g generic grid point 

i inner 

max maximum value 

n nondimensionalized quantity 

o outer boundary 

s surface 

u upstream grid limit 

Procedure 

The geometry input is normally in the wave-drag, 
or Harris, format. (See ref. 2.) This format describes 
the configuration as a set of separate components and 
is convenient for making changes in geometry param- 
eters. The geometry is read into a program (ref. 3) 
that completes the geometry by computing the wing- 
fuselage and fin-fuselage intersection lines and that 
generates a surface grid. (See fig. 1.) To initiate 
the volume grid calculation, a new surface grid is 
computed by interpolating in the initial grid along 
x = Constant lines. Then an outer boundary shape 
is specified, and the grid is generated by interpolation 
between the surface and this outer boundary. Details 
of this procedure follow. 

Surface Grid 

For the initial surface grid lines, points are distrib- 
uted along x = Constant cuts on the fuselage ahead 
of the wing. In the wing region, the wing lofting lines 
are extended to the fuselage surface by the method 
of reference 3 and then are continued around the 
fuselage along x = Constant lines. (See fig. 2.) 

Now, a new surface grid is computed by inter- 
polating in the original grid at x = Constant sta- 
tions. Because the grid is already laid out this way 
on the fuselage, the interpolation there involves only 
the x variable. On the wing, the interpolation is a 
matter of locating where the lofting lines intersect 
the x = Constant planes. 

Care must be taken in distributing the points 
on these lines so that grid continuity is ensured. 
The main problem locations are at the wing and 
fin leading and trailing edges. The grid variable in 
the circumferential direction at x = Constant is de- 
noted 77. (See fig. 3.) To prevent the constant-77 sur- 
face lines from running from the lower to the up- 
per surface or vice versa, the entire configuration is 
treated as separate upper and lower surfaces. Thus, 
each x = Constant line consists of two distinct seg- 
ments, one on the upper surface and the other on the 
lower surface. In the wing region the break is taken 



at the wing tip. On the fuselage ahead of the wing, 
it is taken at the same angular location as the point 
where the wing leading edge intersects the fuselage. 
(See fig. 3.) This system ensures grid continuity in 
the 77 direction at the wing leading edge. 

Grid points are distributed along each of the up- 
per and lower surface x = Constant lines as a func- 
tion of arc length. The arc is initiated at the y — 0 
point, and the cumulative arc length s is computed 
at successive points out to the break point. This 
procedure enables us to express the grid points on 
the line as functions of s (i.e., y(s),z(s)). Then this 
distribution of points can be replaced through inter- 
polation by any desired distribution. At present this 
distribution function is specified by a cubic polyno- 
mial with arbitrary coefficients to permit flexibility 
in designing the surface grid. To ensure continuity in 
the longitudinal direction, this distribution is spec- 
ified globally so that if, for example, the points are 
congested near the wing tips, they will also be con- 
gested on the fuselage ahead of the wing in the region 
that corresponds to the wing tips. An example of the 
resulting surface grid is shown in figure 4. 

At the wing trailing edge a problem arises when 
the x = Constant lines are interpolated from the 
initial surface grid. If the trailing edge is swept, a gap 
occurs between the intersection of an x — Constant 
plane with the fuselage and its intersection with the 
wing. (See fig. 5.) This problem is treated by 
including the gap region as part of the surface grid 
(fig. 6 ) but labeling the grid points in this region 
as wake points. To maintain continuity, the wake is 
assumed to intersect the fuselage at the same angular 
location as the wing trailing edge. A vertical fin 
would be treated similarly if present. 

Volume Grid 

After the surface grid has been established, the in- 
dividual two-dimensional grid surfaces for an Euler- 
marching grid are generated by specifying a corre- 
sponding distribution of points on an appropriate 
outer boundary line and then interpolating between 
the surface points and outer boundary points. An ap- 
propriate outer boundary shape is a circle that lies 
slightly outside the Mach cone whose vertex is at or 
slightly ahead of the nose of the configuration. If 
r s denotes a vector representation of the cut through 
the body surface and r G denotes the outer bound- 
ary shape, then a straightforward transfinite inter- 
polation for the intermediate grid-line shapes (lines 
of constant £) would result from a variation of the 
form 

r g =(1 - a)r s + ar 0 (1) 


where a is a blending, or homotopy, function of £ 
that varies continuously and monotonically from 0 
to 1 . 

For this type of grid the skewness can be severe, 
as is illustrated in figure 7. Somewhat better control 
over the grid shape can be realized by the method of 
reference 4. The inner and outer boundary shapes 
arc nondimensionalized, thereby reducing them to 
the same scale. Then the homotopy function can be 
used to establish a pure shape variation, whereas the 
scaling is specified by an independent function A(£). 
Thus, 

T g = Ap - a) r ns + 0cT no ] (2) 

However, for difficult shapes such as that shown in 
figure 7, this generalization is not sufficient. The ra- 
dial (77 — Constant) grid lines near the wing- fuselage 
juncture need a larger \y\ component to be more 
nearly normal at the body surface. This effect 
could be accomplished by replacing the circular outer 
boundary with a highly eccentric ellipse. However, 
an eccentric ellipse would not be appropriate for flow 
calculation, and it would not help the problem at 
the wing tips, where the 77 — Constant lines require 
a larger z component to reduce the skewness. 

An outer boundary shape that would reduce the 
skewness at both the wing tips and the wing-fuselage 
juncture is a superellipse, which is defined by 

y — xp\ cos 0j 9 sign(cos 6) (3) 

z — p\ sin0) 9 sign(sin 0)+ z (4) 

where \ I s the eccentricity parameter and q is the 
superellipticity exponent. This outer boundary 
shape does reduce the skewness at the body sur- 
face, but it does not yield a practical grid away from 
the surface. However, this problem can be overcome 
by making a further generalization of the transfinite 
interpolation formula. Returning to equation ( 1 ), we 
replace the fixed circular outer boundary shape r 0 
with a function that varies parametrically with £ such 
that it becomes a circle at the largest value of £. At 
small values of £ (near the surface), t g becomes the 
superellipse of equation (3). Between these two ex- 
tremes, t g is defined as a smooth blend of the inner 
function r l0 and the outer function by a homotopy 
function 7 (£): 

r o =(1 - l)Tio + 1^00 (5) 

This function replaces r G in equation (1). 
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In equation (3), both the eccentricity param- 
eter x and superellipticity exponent q are speci- 
fied in terms of the cross-section aspect ratio 
A = ymax/Fuselage semiheight. Thus, 

x = l + p & t) <6) 

where p is an adjustable input parameter of the order 
of 1 and 

q( C)= V>+(1 -#7(0 ( 7 ) 

where the variable 

v- = 4 ( 8 ) 

Also in equation (4) the z symmetry line z = z is 
taken where the superellipse is intersected by a line 
drawn from the fuselage center through the wing tip. 

Thus, near the surface, the influence of r 0 on the 
grid shape is that of a superellipse, but far from the 
surface it becomes that of a circle. This effect is seen 
in the example shown in figure 8. 

Several typical surfaces of the complete three- 
dimensional grid are shown in figure 9. Figure 10 
shows results of an Euler flow calculation on this grid. 
The marching procedure used for this calculation is 
similar to that of reference 5. For the half plane 
y > 0, the grid dimensions are 97, 97, and 49 in the 
f , Tj, and C directions, respectively. 

Subsonic Grid 

A grid for calculating subsonic flows over a con- 
figuration must extend well ahead of, behind, and 
to the sides of the configuration. One type of sub- 
sonic grid can be generated starting with the same 
techniques used to generate the Euler grid but by 
changing the outer boundary shape from a cone to a 
more appropriate shape. Then, the two-dimensional 
surfaces at x = Constant cuts are changed into three- 
dimensional surfaces by distorting them out of the 
x = Constant planes. 

An example has been calculated for which the 
outer boundary is taken to be a semiellipsoid up- 
stream that is connected to a cylinder downstream 
at the x station for which the wing attains its max- 
imum span. (See fig. 9(a).) The radius distribution 
of the ellipsoid is 


,, / 1 

'(x - x s ) ' 

2 

r max W 1 ~ 

_(x u — x s )_ 



where x u is the grid station that is farthest upstream. 
The radius of the downstream cylinder is, of course, 

T — ^max ( 10 ) 

The two-dimensional x = Constant grid surfaces are 
converted to three-dimensional surfaces by supple- 
menting equation (2) with an x variation 

x g =(1 - 6)x s -I- Sx 0 (11) 

where the subscript g is the generic grid point and <5 
is a homotopy function. Typical surfaces of this grid 
are shown in figure 11. 

Concluding Remarks 

Several procedures for the automatic generation 
of flow computation grids for relatively simple config- 
urations have been described. For supersonic flows, a 
quasi-two-dimensional grid for marching Euler codes 
was developed, and some sample results in graphi- 
cal form were included. For subsonic flows, the pro- 
cedure was modified by distorting the x — Constant 
grid surfaces out of plane and by specifying a more 
appropriate outer boundary shape. The techniques 
are algebraic and are based on a generalization of the 
transfinite interpolation method. 


NASA Langley Research Center 
Hampton, VA 23681-0001 
May 16, 1994 

References 

1. Smith, Robert E., ed.: Software Systems for Surface 

Modeling and Grid Generation. NASA CP-3143, 1992. 

2. Middleton, W. D.; Lundry, J. L.; and Coleman, R. G.: 
System for Aerodynamic Design and Analysis of Super- 
sonic Aircraft . Part 2- User’s Manual. NASA CR-3352, 
1980. 

3. Barger, Raymond L.; and Adams, Mary S.: Auto- 

matic Computation of Wing-Fuselage Intersection Lines 
and Fillet Inserts With Fixed- Area Constraint NASA 
TM-4406, 1993. 

4. Barger, Raymond L.; and Adams, Mary S.: Semianalytic 
Modeling of Aerodynamic Shapes. NASA TP-2413, 1985. 

5. Morrison, J. H.; and Korte, J. J.: Implementation of 
Vigneron’s Streamwise Pressure Gradient Approximation 
in Parabolized Navier-Stokes Equations. AIAA J., vol. 30, 
no. 11, Nov. 1992, pp. 2774-2776. 


3 




(b) Configuration with computed intersection line (lower surface). 
Figure 1. Calculated wing- fuselage intersection. 
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Figure 2. Close-up view of intersection line and surface grid near wing leading edge. 
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(a) Grid section ahead of wing leading edge intersection with fuselage. 



(b) Grid section behind wing leading edge intersection with fuselage. 

Figure 3. Grid sections displaying grid coordinates and circumferential distribution of grid lines. 


Figure 4. Revised surface grid to initialize Euler-marching grid. x = Constant lines are computed by 
interpolating in original grid and redistributing points. 



Figure 5. Plan view illustrating gap that occurs in x = Constant grid line when trailing edge is swept. 
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(a) Wing-body juncture region. 


Figure 7. Grid skewness problem at body surface. 
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Figure 8. Grid skewness reduction at body surface. 
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(b) Wing-tip region. 
Figure 8. Concluded. 
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(a) j = 1 and 97 surfaces. 


Figure 11. Typical surfaces of subsonic flow grid. 
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(b) j = 49 surface. 
Figure 11. Continued. 



(c) k = 12 surface (grid diluted to improve clarity). 
Figure 11. Continued. 
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(d) i = 24 surface. 
Figure 11. Concluded. 
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